Quantifying Wrinkle Features of Thin Membrane Structures 


Mindy B. Jacobson” 

NASA Goddard Space Flight Center, Greenbelt, Maryland, 20771, USA 

Takashi Iwasa + 

JAXA, Institute of Space and Astronauticai Science, Sagamihara, Kanagawa, Japan 

and 

M.C. Natorr 

JAXA, Institute of Space and Astronauticai Science, Sagamihara, Kanagawa, Japan 


For future micro-systems utilizing membrane based structures, quantified predictions of 
wrinkling behavior in terms of amplitude, angle and wavelength are needed to optimize the 
efficiency and integrity of such structures, as well as their associated control systems. For 
numerical analyses performed in the past, limitations on the accuracy of membrane 
distortion simulations have often been related to the assumptions made. This work 
demonstrates that critical assumptions include: effects of gravity, supposed initial or 
boundary conditions, and the type of element used to model the membrane. In this work, a 
0.2m x 0.2m membrane is treated as a structural material with non-negligible bending 
stiffness. Finite element modeling is used to simulate wrinkling behavior due to a constant 
applied in-plane shear load. Membrane thickness, gravity effects, and initial imperfections 
with respect to flatness were varied in numerous nonlinear analysis cases. Significant 
findings include notable variations in wrinkle modes for thickness in the range of 50 pm to 
1000 pm, which also depend on the presence of an applied gravity field. However, it is 
revealed that relationships between overall strain energy density and thickness for cases 
with differing initial conditions are independent of assumed initial conditions. In addition, 
analysis results indicate that the relationship between wrinkle amplitude scale (W/t) and 
structural scale (L/t) is independent of the nonlinear relationship between thickness and 
stiffness. 
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I. Introduction 

W ITH the need to develop increasingly more efficient and lighter mechanical and structural systems for 
spacecraft, many researchers have for decades been experimenting with various types of polymer membranes 
and working to simulate their behavior. 1 " To date, the most detailed simulations have been limited to focusing on 
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macro-scale global distortion and the overall stress field: analysis tools have not been robust enough to quantify 
wnnkie features such as amplitude, wavelength and orientation angle. In the future, quantified predictions of 
membrane wrinkling behavior in terms of these three characteristics can be used to optimize the efficiency and 
integrity of membrane structures and associated control systems. 

Historically prominent space-based examples of structural systems employing thin membrane materials include 
the Echo-1 passive communications satellite, high altitude balloons, and Fresnel lens reflectors. Unlike these 
systems, future membrane based space structures are likely to 
require controllability for various precision applications. Because 
controllability of wrinkled surfaces is a highly complex problem, the 
current design concept for such structures enforces precision 
controls on the fabrication processes, rather than on real-time 
functionality. Examples of some current design concepts of future 
space-based optics include the James Webb Space Telescope 
(JWST) sun shield, Constellation-X (CON-X) mirrors. Extreme 
Universe Space Observatory (EUSO) lenses, and ultra-light 
precision optics, as shown in Figure 1 (a)-(d), respectively. Each of 
these structures is highly sensitive to wrinkle patterns in only one 
dimension, which is an order of magnitude smaller than the other 
two planar dimensions. Controllability of this dimension and 
distortion in the out-of-plane direction significantly impacts the 
function of each system. 

Lessons learned from each of these classes of membrane 
structures can be applied to future controllers and micro-systems 
only as quantifiable analysis tools become available. The analytical 
study documented herein includes wrinkle simulations for the basic 
case of a square membrane fixed at each comer and subjected to in- 
plane shear. Use of the Mixed Interpolation of Tensorial 
Components (M1TC) shell element has enabled this detailed 
quantitative study for reasons discussed in the following section. The subsequent section describes numerous 
nonlinear analysis cases considered, in which membrane thickness, gravity effects, and initial imperfections with 
respect to flatness were varied. Simulation results are presented and quantified in terms of wrinkle amplitude, 
orientation and wavelength. 

Results of this work are significant as a tool for: ( 1 ) assessing numerical analysis sensitivities in preparation for 
data correlation, and (2) developing alternative tools for analyzing membrane wrinkling behavior predicted by 
experiments or finite element methods. ' 



Figure 1. Space-based membrane 
structures: (a) JWST suits hie Id 
concept, (b) CON-X mirror, (c) EUSO 
lenses, (4) ultra-light precision optics 


II. Wrinkling Theory' 

Although a complete history of the development of theories for explaining wrinkling behavior are beyond the 
scope of this paper, three theories traditionally used to study and simulate wrinkled membranes are briefly 
summarized here. 

First, Tension Field Theory (TFT) mathematically describes wrinkling in terms of solid mechanics. 4 It asserts 
that membranes can be described as structures of sufficiently little thickness such that they have zero bending 
stiffness. As a result, they have no capacity to support compressive load or carry compressive stress. 

Plate Theory is assumed when considering structural responses that are dominated by in-plane structural 
stiffness. Out-of-plane, or transverse loads are considered negligible by comparison. Elements treated according to 
this theory are not expected to undergo significant amounts of bending, although they do have significant thickness. 

Classic Shell Theory is assumed when a thin structural element is loaded both in-plane and out-of-plane, and 
bending is expected. The theory includes non-negligible bending stiffness and can be used for considering nonlinear 
deformation, or strain. Further, it considers that the shell is an inextensible surface and is modeled at the mid-plane 
of that surface - therefore, transverse stresses through the thickness of a surface are not captured. 

For elastic materials, such as membranes considered in this work, inextensibility is a fundamental principle 
stating that the length of a distorted path cannot exceed the length of an undistorted path of the same structure. In 
other words, geometric compatibility between initial and final states must be preserved. For wrinkling phenomena, 
this principle relates to the transverse direction. Hence, Plate Theory is insufficient for characterizing wrinkling 
behavior. Classic Shell Theory is sufficient only when formulations including transverse shear stresses are 
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employed. Such formulations include, for example. Discrete Kirchhoff Shells and Mixed Interpolation of Tensoria! 
Components (M1TC) shells. 

The MITC sheii element was developed by K. j. Bai’ne and E. N. Dvorkin to combine the plate and shell i’neories 
in finite element analyses . 5 ' h ' 1 The stress and strain tensors at each element node are interpolated along each edge of 
the element, as well as through the element thickness. The interpolation scheme used makes this element more 
resistant to shear locking than other shell elements. Altogether, the characteristics of the MITC shell element allow 
them to capture true strain energy density levels. This is critical, given that 
wrinkling is best described as the minimum energy response for 
accommodating in-plane strain while also satisfying the inextensibility 
principle in the transverse direction. 

Mathematically, wrinkles have been described in terms of both TFT 
and shell theory. For the purposes of this study, to be consistent with 
MITC element formulation, the principles of shell theory are used to 
characterize wrinkle states, as described in the following Table 1. 


Table!: Wrinkle State Criteria 
Defined by Principal Stresses 


Taut 

Wrinkled 

Slack 


On > 0 and O 21 > 0 
an > 0 and 022 ^ 0 
an < 0 and 022 ^ 0 


III. Analysis 


A. Method 

The computational tool used to perform these simulations is a Finite Element Analysis Program: Personal 
Version (FEAPpv). 1 ' 9 11 The 4-node MITC shell element was added by authors to the element library'. The Newton 
and modified Newton methods are used to compute nonlinear analysis solutions. 

B. Model 

The 0.2m x 0.2m membrane finite element model is 
showti in Figure 2. The mesh was discretized with regard 
to optimal aspect ratio for best performance of the 4-node 
MITC shell element. 11 Computation time is of secondary 
importance in this study; the model consists of 10000 
elements (0.002rrf each) and 10201 nodes. 

The membrane is constrained at 3 of 4 comers in all 
degrees of freedom except out of plane rotation. At the 
fourth comer of the membrane, a small displacement is 
applied in the direction of the +X axis (8x = 0.001m) and 
constraints are applied against translation in the Y and Z 
directions, as w'eli as against rotation about the X and Y 
axes. 

The magnitude of applied displacement was set with 
consideration for the assumption that this analytical study 
is based on linear elastic behavior of the membrane. Resultant global stress levels are not to exceed the 5% yield 
point for H-type Kapton; refer to Table 2 for assumed material properties. 


Table 2: H-type Kapton Material Properties 


Property 

Value 

Unit 

Tensile Strength, Ultimate 

243.04 

MPa 

Yield Point (3% offset) 

95.06 

MPa 

Yield Point (5% offset) 

123.48 

MPa 

Elongation, Ultimate 

80 

% 

Young’s Modulus, E 

2959.6 

MPa 

Poisson's Ratio, v 

0.34 


Density, p 

1.42 

gem’ 



Figure 2: Finite Element Model 
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C. Cases 

Parameters discussed in the previous section are held 
constant for all numerical simulations. Each analysis case 
performed in this study was derived based on the need to 
systematically evaluate effects of additional analytical 
variables. Table 3 lists the variable parameters treated. 

For each thickness, geometrically nonlinear analyses 
were performed in order to assess the winkling behavior 
both with and without the gravity effect, which is assumed 
perpendicular to the plane of the membrane. For the cases 
including this effect, membrane sag due to gravity is taken 
to be the initial distorted shape prior to applying the shear 
load (as an enforced displacement, 5x = 0.001 m). 

For cases not including the gravity effect, the initial flatness of the membrane is assumed to not be perfect. That 
is, in reality, no square membrane surface can be considered perfectly flat when it is suspended at each of its 4 
comers. In this study, the initial deformed shape (or, flatness imperfection) is determined in one of two ways. The 
first method is to apply random values of out-of-plane distortion at each node. These values are ±10% of the 
membrane thickness and are taken from a normal random number (NRN) generator, and they are applied as initial 


Tabic 3l PaTauicici 5 

Varied for Analytical Study 

Parameter 

Values 

Units 

Thickness, t 

50, 125, 150, 175, 
300, 500, 1000 

pm 

Gravity 

0 or (p-At)-G 

M-N 

Initial Flatness 
Imperfection 

- 1 0%t < NRN < 
+10%t 

pm 


Table 4: Description of Analysis Cases 


Case 

Applied 

Load 

Initial Condition 

Boundary 

Condition 

Illustration 

; “Self-Weight” 

Effective 
Gravity 
Force at 
Each Node 

None 

5DOF Constrained 
at All 4 comers 

i; 

4 » 

„x 

DOS 

0 Q j 

“Gravity” 

8x=0.001m 
at one 
comer 
(X=0.2m, 
Y=0.0m) 

Restart of Self- 
Weight Case 

5DOF Constrained 
at 3/4 comers; 
4DOF constrained 
at comer of 
applied load 

13 

•j ■. 

“No Gravity 
NRN” 

8x=0.00 1 m 
at one 
comer 
(X=0.2m, 
Y=0.0m) 

Flatness 

Imperfection due 
to random 
distribution of 
distortions (8z) 

5DOF Constrained 
at 3/4 comers; 
4DOF constrained 
at comer of 
applied load 

-V' . 

**■*•*•» ei - 

\ / 

■■ DK 

“No Gravity 
S-W” 

Sag shape due to 
Self-Weight 
analysis 

1 

I 

“J\ ~ ..Vs — - 

0*' 
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conditions. Alternatively, the initial deformation is assumed to correspond to the sag shape due to gravity, but 
applied as a geometric initial condition, rather than taken as the result of an applied load. 

Table 4 summarizes each of the analysis cases considered in this study. All cases, but the first ("Self- Weight”), 
are considered to be wrinkling analysis cases for which a detailed study of the results has been conducted. The 
"Self-Weight” case is taken to be a baseline for comparing the quantified effects due to gravity simulated in the 
other three wrinkling cases. Note that each membrane thickness given in Table 3 is considered for all cases. 

IV. Results 


A. Self-Weight Case 

The purpose of this case is to simulate the static 
effect of a (1G) gravity field on a square membrane 
held fixed (5DOF) at all four comers. The gravity field 
is simulated by computing an effective gravity force 
(F g ) for each membrane thickness, which is evenly 
distributed over the entire area (A) and applied 
discretely at each node. 

Figure 3 illustrates the resultant deformed shape 
that is characteristic for each membrane thickness, 
from 50 pm to 1000 pm. Note that the shape of any 
cross section profile is not a catenary. 

For comparison with the three wrinkling analysis 
cases, total strain energy density (U) for each 
membrane thickness is presented in Figure 4; a non- 
dimensional structural scale (L/t) parameter is given 
along the abscissa. 





Figure 3: Characteristic Deformation Due to Gravity 


It is clear that U rapidly decreases with thickness. The rate of decrease is related to the rate of decrease in 
membrane stiffness, which is a compound function of both t and t 3 for in-plane and bending stiffness, respectively. 
Overall work required to distort the membrane is significantly less for small thickness than for large: therefore, this 
relationship is expected. 
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Figure 4: Strain Energy Density Due to Gravity Sag 


B. Gravity Case 

The solution of this case is obtained from a nonlinear solution re-start of the “Self- Weight” case while applying 
the comer shear force. In this sense, the simulation may be considered to be the result of an initial applied load, or 
force. (This is contrary to the “No Gravity S-W” case, described below, in which the gravity effect is modeled as a 
geometric initial condition, or distortion.) 

Although simulations for each thickness given in Table 3 were performed, for brevity only the deformed shapes 
for t = 50, 175, 300, and 1000 pm are illustrated in Figure 5. Simulation results indicate that only when t < 500um, 
maximum out-of-plane distortions (8z, in either the +Z or -Z direction) decrease as thickness decreases. However, 
there is no clear thickness dependent trend for maximum stress levels, which occur at the point of application of the 
comer shear load; the difference between maximum and minimum values is only 5.7%. 
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Element-by-element analysis of wrinkle state, as 
defined in Table 1, for each membrane thickness has 
revealed that for t > 300jj.ni, wrinkle cniena based on 
major principal stresses does not accurately 
characterize membrane distortion because some 
regions do carry some level of compressive stress 
without buckling due to sufficiently high levels of 
bending stiffness. This may be confirmed by 
examining the deformation modes given in Figure 5; 
true wrinkling behavior (a result of localized buckling) 
is dominant only when t < 300jim. 

C. No Gravity NRN Case 

This analysis case was prepared by enforcing, as a 
geometric initial condition, small out-of-plane (5z) 
distortions at each node location. The magnitude of 
these distortions is not greater than 10% of the 
membrane thickness, and they were taken from a 
normal random number (NRN) generator. 

This initial condition was chosen in order to induce 


t = 50 pm 
Max. jSzj - 2.4 mm 
CTpmax= 141.1 MPa 
Max. |x X y| = 55.7 MPa 

! 

|KjKPH| I 

t= 175 pm 
Max. Szj = 2.8 mm 
Optox = 142.8 MPa 
Max. jXjjyj = 56.9 MPa 

| 


t - 300 pm 
Max. |5z| = 2.8 mm 
Opmax = 144.7 MPa 
Max. |x X y| = 57.6 MPa 


t = 1 000 pm 
Max. |5zi = 2.6 mm 
Opmax = 143.8 MPa 
Max. !x xy | = 57.8 MPa 

J! 

8 0 



Figure 6: Wrinkled Membranes Due to Applied 
Corner Shear Load, Including Initial NRN 
Imperfection 


t = 50 pm 

Vlov IXt! ~1 mm 

Opmax = 141.4 MPa 
Max. jx xy | = 55.8 MPa 

! 


t = 1 75 pm 
Max. |5z| = 2.9 mm 
Opmax = 142.8 MPa 
Max. |x xy | = 56.9 MPa 

i”: 


t = 300 pm 
Max. |8z| = 3.0 mm 
Opmax = 144.5 MPa 
Max. |x xy | = 57.5 MPa 

a 

2 c. 


t - 1 000 pm 
Max. |5z = 2.6 mm 
Opmax = 1 50.0 MPa 
Max. |X xy ! = 59. 1 MPa 

i": 

•* 

** — IH 0,4 


Figure 5: Wrinkled Membranes Due to Applied 
Corner Shear Load, Including Gravity Field 


a perturbation in the flatness of the membrane - in the 
absence of a gravity effect. The perturbation is 
necessary because the assumption of a perfectly flat 
surface is an inherently stable system. The nature of the 
NRN displacement initial condition allows for the 
computation of a solution based on imperfect 
bifurcation, meaning that no single eigenvector 
trajectory is followed through the solution sequence. 

As in Figure 5 for the "Gravity’' case. Figure 6 
illustrates wrinkled membranes only for t = 50, 175, 
300. and 1000 |im. A quantitative comparison of the 
stresses and deformations associated with these results 
and those given in Figure 5 reveals no significant 
differences. However, there are obvious qualitative 
differences as thickness decreases. 

Therefore, in order to truly characterize wrinkling 
behavior it is necessary to quantify wrinkle orientation 
in addition to amplitude. This has been done, for this 
study, in a rudimentary way by plotting maximum and 
minimum principle stresses at a constant cross section 
in the membrane, for each thickness' as shown in 
Figure 7. The most significant detail to notice from this 


This topic is expanded and thoroughly treated in a separate study entitled “A Spectrum Evaluation Method for 
Wrinkled Membranes Using Two-Dimensional Fast Fourier Transform” by T. Iwasa, et. al. 3 


6 

American Institute of Aeronautics and Astronautics 









comparison is that minor principle stresses tend to reverse, and with greater frequency, as thickness decreases. This 



Figure 7: On and 022 Along Membrane Cross Section at Y = 0.14m, including Initial NRN Imperfection 


is the cause of increasingly numerous and skewed wrinkles when thickness is less than 300 pm. 

D. No Gravity S-W Case 

This wTinkling case differs from those presented in 
the previous two sections in that the assumed initial 
condition is determined by the appropriate degree of 
gravity sag deformation, depending on thickness, that 
was computed in the “Self- Weight” case. This is 
unlike the “Gravity” analysis case because the initial 
condition is applied as a geometric constraint, rather 
than as a distributed applied load (Fg). 

Figure 8 gives analysis results similar to those 
reported in Figures 5 and 6, for t = 50, 1 75, 300, and 
1000 pm. These results reveal similar trends with 
respect to wrinkle state; when t > 300 pm deformation 
may not be accurately characterized by the term 
"wrinkled”, according to the state definitions given in 
Table 1. In addition, stress levels and out-of-plane 
deformation are not significantly different than the 
other two wrinkling cases. 

However, the deformed state when t = 50 pm is 
different than the previous two cases. In addition, 
differences are discovered when carefully analyzing 
membrane deformation at constant cross sections using 
the wrinkle scale parameter (W/t). This is discussed in 
greater detail in the following section. 

V. Discussion 

In order to visualize how wrinkling behavior is 
effected by the different initial conditions it is 
beneficial to examine cross section plots at constant 
thickness for each analysis case. The non-dimensional 
wrinkle scale parameter W/t is plotted versus structural 
scale L/t, where W=peak out-of-plane displacement 
and L=edge length. For illustration, here, the cross section plots for t = 50, 1 75, 300, and 1000 pm at y = 0.14m are 
shown in Figure 9. 

A detailed examination of such plots for each membrane thickness reveals that deformation mode is independent 
of initial conditions when t < 150pm. Also, for any thickness, the deformation modes for "Gravity” and “No Gravity 


t = 50 pm 
Max. |5zj = 2.7 mm 
Opmax = 138.5 MPa 
Max. ]T xy | = 54.7 MPa 


t= 175 pm 
Max. 8zj = 2.5 mm 
Opmax = 142.9 MPa 

Max. iT xv i = 56.5 MPa 

' 

*Xmrnrnm ”* 

t = 300 pm 
Max. 8zj = 2.5 mm 
Opmax - 145.4 MPa 
Max. T xy j = 57.5 MPa 

•Xmmm •** 

t = 1 000 pm 
Max. 5z| = 3.1 mm 
Ctpmax = 153.3 MPa 
Max. |T xy | = 60.8 MPa 

**■*« 


Figure 8: Wrinkled Membranes Due to Applied 
Corner Shear Load, Including Initial Self-Weight 
Gravity Sag Shape 


7 

American Institute of Aeronautics and Astronautics 












Figure 9: (W/t) Along Membrane Cross Section at Y = 0.14m 


S-W” are the same, but offset by the initial 8z offset due to the gravity sag shape. The disparity between these two 
deformation modes and that for the “No Gravity NRN” case becomes noticeable when 175pm < t < 1000pm and a 
phase shift between wrinkle peaks becomes apparent. But. when t = 1000pm. deformation mode is again 
independent of initial conditions. Another prominent trend to notice is that wavelength increases with thickness, 
regardless of assumed initial conditions. 

Figure 10 illustrates the relationship between the 
wrinkle (W/t) and structural (L/t) scale parameters. 

(Recall that wrinkle scale is defined by the peak out-of- 
plane displacement.) It is important to note that this 
relationship is linear for the “Gravity" case when 
t < 150pm. This indicates that, for certain structural 
scales, analysis results may be scaled for correlation with 
ground-based experimental data - but only for square 
geometries. That is, previous research studying similar 
behavior but for a circular membrane has revealed that the 
same relationship is not linear. 12 Rather, it increases 
nonlinearly towards a constant value as L/t increases. So, 
the relationship between wrinkle scale and structural scale 
depends on geometry, not initial conditions. 

Plotting total strain energy density versus structural 
scale reveals the same relationship for all three wrinkling 
analysis cases as was revealed for the “Self-Weight” case. 
Refer to Figure 1 1 , which shows that volumetric strain 
energy decreases rapidly with thickness and approaches a 
constant value. Therefore, the minimum energy response 
for accommodating in-plane strain (in the form of 
wrinkles) approaches a constant as thickness becomes 
small. However, as previously noted, when t > 500pm 
wrinkling is not an accurate description of the deformation 
due to in-plane shear. 

With respect to strain energy density, it is also 
interesting to note that the ratio of results for the “Gravity" 

Figure 11: Strain Energy Density versus (L/t) 
for AH Wrinkling Cases 
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Figure 10: Overall Maximum (W/t) versus (L/t) 












wrinkling case and the baseline “Self-Weight” cases 
increases rapidly and approaches a constant value as 
thickness decreases (refeT to Figure 12). 

Additional analyses have revealed that the ratio of 
strain energy density in the "Gravity” case to that for 
either of the other two wrinkling cases is 1.0 ± 5%, 
regardless of thickness. Therefore, initial conditions do 
not significantly impact the energy state of a wrinkled 
membrane. Other parameters relating wrinkle amplitude 
to orientation are necessary to fully quantify and qualify 
wrinkle behavior. 
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Figure 12: Strain Energy Density Ratio: 
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VI. Summary 

A systematic analytical study of the affects of thickness, gravity, and initial conditions on membrane wrinkling 
behavior has been presented. The universally applied static load in this study is in-plane shear, applied at one comer 
of a square membrane held fixed at the remaining three comers. This study is unlike others previously conducted on 
wrinkling behavior because of the geometry and numerical methods used to simulate behavior. 

Detailed analysis of the simulation results for three different wrinkling cases reveal important relationships 
between deformation and thickness, as well as initial conditions: 

• As thickness decreases from the maximum considered in this study (1000 pm), the global wrinkling 
behavior is sensitive to initial conditions only when thickness is less than 1 50p.m. 

• As thickness increases, greater than 300 pm, membrane deformation may no longer be mathematically 
characterized as wrinkling. 

• Overall stress and strain levels are not significantly effected by assumed initial conditions, though wrinkle 
orientation is. 

• Wrinkle wavelength increases with thickness for any assumed initial condition. 

• By plotting wrinkle scale with respect to structural scale, it is shown that wrinkle scale is independent of the 
nonlinear relationship between thickness and stiffness. 

• The relationship between wrinkle scale and structural scale is linear when t < 1 50 pm for the “Gravity” case 
in this study. 

• The ratio of overall strain energy density in the "Gravity” case to that for either of the other two wrinkling 
cases is approximately 1 .0. 

All of these results may be useful for correlating test data from a ground-based experiment. In addition, they 
reveal the fact that quantifying wrinkle behavior due to in-plane shear is a laborious and time-consuming task when 
utilizing finite element simulations. Because knowledge of precise wrinkle characteristics such as orientation and 
wavelength may be necessary for optical applications of membrane structures, this analytical study may be useful as 
the background for developing alternative wrinkling analysis tools capable of predicting such details outright. 
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